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ABSTRACT 

We study the dynamics of stellar-mass black holes (BH) in star clusters with particu- 
lar attention to the formation of BH-BH binaries, which are interesting as sources of 
, gravitational waves (GW). In the present study, we examine the properties of these 

BH-BH binaries through direct N-body simulations of star clusters using the NBODY6 
code on Graphical Processing Unit (GPU) platforms. We perform simulations for star 
clusters with $C 10 5 low-mass stars starting from Plummer models with an initial 
O ' population of BHs, varying the cluster-mass and BH-retention fraction. Additionally 

we do several calculations of star clusters confined within a reflective boundary mim- 
icking only the core of a massive star cluster which can be performed much faster 
than the corresponding full cluster integration. We find that stellar-mass BHs with 
masses ~ 10M® segregate rapidly (~ 100 Myr timescale) into the cluster core and 
form a dense sub-cluster of BHs within typically 0.2 — 0.5 pc radius. In such a sub- 
cluster, BH-BH binaries can be formed through 3-body encounters, the rate of which 
can become substantial in dense enough BH-cores. While most BH binaries are finally 
ejected from the cluster by recoils received during super-elastic encounters with the 
' single BHs, few of them harden sufficiently so that they can merge via GW emission 

within the cluster. We find that for clusters with N > 5 x 10 4 , typically 1 - 2 BH-BH 
mergers occur per cluster within the first ~ 4 Gyr of cluster evolution. Also for each 
of these clusters, there are a few escaping BH binaries that can merge within a Hubble 
time, most of the merger times being within a few Gyr. These results indicate that 
intermediate-age massive clusters constitute the most important class of candidates 
for producing dynamical BH-BH mergers. Old globular clusters cannot contribute sig- 
nificantly to the present-day BH-BH merger rate since most of the mergers from them 
would have occurred much earlier. On the other hand, young massive clusters with 
ages less that 50 Myr are too young to produce significant number of BH-BH mergers. 
We finally discuss the detection rate of BH-BH inspirals by the "LIGO" and "Ad- 
vanced LIGO" GW detectors. Our results indicate that dynamical BH-BH binaries 
constitute the dominant channel for BH-BH merger detection. 

Key words: gravitational waves - stellar dynamics - scattering - methods: N-body 
simulations - galaxies: star clusters - black hole physics 
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1 INTRODUCTION 

Star clusters, e.g., globular clusters (henceforth GC), young 
and intermediate-age massive clusters and open clusters har- 
bor a large overdensity of compact stellar remnants com- 
pared to that in the field by virtue of their high density and 
the mass segr egation of the compac t remnants towards the 
cluster core (|Hut fc Verbuntl Il983t ). These compact stars, 
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which are neutron stars (henceforth NS) and black holes 
(henceforth BH), are produced by the stellar evolution of 
the most massive stars within the first ~ 50 Myr after clus- 
ter formation. These compact stars, being generally heav- 
ier than the remaining low-mass stars of the cluster, segre- 
gate quickly (within one or a few half-mass relaxation times) 
to the cluster core, forming a dense sub-cluster of compact 
stars. Such a compact-star sub-cluster is of broad interest 
as it efficiently prod uces compact-star binaries through dy - 
namical encounters l|Hut fe Verbunj|l983l : [Hut etaU 1992), 
e.g., X-ray binaries and NS-NS and BH-BH binaries, which 
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are of interest for a wide range of physical phenomena. For 
example, X-ray binaries are the primary sources of GC X-ray 
flux, while mergers of tight NS-NS binaries (through emis- 
sion of gravitational waves) is the most likely scenario for 
the production of short duration GRBs and both NS-NS and 
BH-BH mergers are very important sources of gravitational 
waves (henceforth GW) detectable by future gravitational 
wave ob servatories like "Advanced LIGO" (AdLIGO) and 
"LISA" (|Amaro-Seoane et ai] |2007). Double-NS systems are 
also interesting because they c ould be observable as double- 
pulsar systems llRansoml [2 0081, allowing important tests of 
General Relativity (|Kramerll200a ). In the present work, we 
investigate the dynamics of stellar-mass BHs in star clus- 
ters, with particular emphasis on dynamically formed BH- 
BH binaries. Such binaries are strong sources of GWs as 
they spiral-in through GW radiation, a process detectable 
out to several thousand Mpc distances. 

Tight BH binaries that can merge within a Hub- 
ble time can also be produced in the galactic field from 
tight stellar binaries, which r esult from stellar evolu- 
tion of the components ( e.g., [B ulik & Bclczynski 120031 ; 
lO'Shaughnessv et al.l 120071 . |200S| ; iBelczvnski et al.1 l2007h . 
However, Belczvnski et all (|2007t) have shown with revised 



binary-evolution models that the majority of potential BH- 
BH binary progenitors actually merge because of common 
envelope (CE) evolution which occurs when any of the bi- 
nary members crosses the Hertzsprung gap. They found that 
this reduces the merger rate by a factor as large as ~ 500 
and the resulting AdLIGO detection rate of BH-BH binary 
mergers in the Universe from primordial binaries is only ~ 2 
yr _1 , subject to the uncertainties of the CE evolution model 
that has been incorporated. Note that the corresponding NS- 
NS merger rate, as the above authors predict, is considerably 
higher, ~ 20 yr _1 . In view of the above result, the majority 
of merging BH binaries are possibly those that are formed 
dynamically in star clusters. 

As studied by seve ral authors earlier (|Merritt et al.l 
l200fj ; lMackev et al.ll2007u , black holes, formed through stel- 
lar evolution, segregate into the cluster core within ~ 0.3 pc 
and form a sub-cluster of BHs, where the density of BHs is 
large enough that BH-BH binar y formation through 3-body 
encounters becomes important (|Heggie fc Hutll2003l ). These 
dynamically formed BH binaries then "harden" through re- 
peated s_u£er-ejBstic_«rTOUirte£S w^Hi surrounding BHs 
(|Heggielll975l ; lBaneriee fc Ghoshll2006l ). The binding energy 
of the BH binaries released is carried away by the single and 
binary BHs involved in the encounters. This causes the BHs 
and the BH binaries to get ejected from the BH-core to larger 
radii of the cluster and they heat the cluster while sinking 
back to the core through dynamical friction l|Mackev et al.l 
120071 ). Most of the energy of the sinking BHs is deposited in 
the cluster core as the stellar density is much higher there. 
As the BH binaries harden, the encounter-driven recoil be- 
comes stronger and finally the recoil is large enough that 
the encountering single BH and /o r the BH binary e scape 
from the cluster (see Sec. |3j also iBenaccmistal [2002) . Be- 
cause of the associated mass-loss from the cluster core, this 
also results in cluster heating. These heating mechanisms 
result in an expansion of the cluster , as st u died in detail 
by sev e ral authors, e.g . , Merritt et al.l (|2006l ) ; iMackev et al.l 
<|2007h . IMackev et all ([20071 ) found notable agreement be- 
tween the core expansion as obtained from their N-body 



simulations with the observed age-core radius correlation 
for star clusters in the Magellanic Clouds. 

The dynamics of stellar mass BHs and formation 
of BH binaries in star clusters and the resulting rate 
of GW emission from close enough BH binaries has 
been studied by several authors us ing N-body integrations 
jPortegies Zwart fc McMillanll2000l) . nume rical 3-body scat- 
tering ex periments ( Giiltckin ct.al. 2004) or Monte-Carlo 
metho ds (jO'Learv et.a]||2006h . |Portegies Zwart fc McMillanl 
(2000) considered the rate of mergers of escaping BH bina- 
ries from various stellar systems, e.g., massive GCs, young 
populous clusters and galactic nuclei. They estimated the 
merger rates within 15 Gyr (their adopted age of the Uni- 
verse) by assuming the binding-energy (Et) distribution of 
the escaping BH binaries to be uniformly distributed in 
log Eb, as inferre d from simulations of N ~ 2000 o r 4000 
star clusters (see iPortegies Zwart fc McMillanl l2000l for de- 
tails). Considering the space-densities of the different kinds 
of star clusters, they found that while for LIGO the cor- 
responding total (i.e., contribution from all types of star 
clusters) detection rate i s negligible, it c an be as high as 
~ 1 day" 1 for AdLIGO. iGiiltekin etlail (|2004) performed 
sequential numerical integrations of BH binary-single BH 
close encounters in a uniform stellar background, where, in 
between successive encounters, the BH-binaries were evolved 
due to GW emission. From such simulations, these authors 
studied the growth of BHs through successive BH-binary 
mergers for the first time. In a more self-consistent study 
of the growth of BHs and the possibility of formation of in- 
termedia te mass black h oles t hrough successive BH-binary 
mergers, lO'Learv et.all (|2006T ) used the Monte-Carlo ap- 
proach and considered "pure" BH clusters that are dynam- 
ically detached from their parent clusters which ca n be ex- 
pecte d to form due to mass stratification instability l|Spitzerl 
ll987l V These authors considered dynamically formed BH 
binaries from 3-body encounters in such BH clusters uti- 
lizing thfioretic^Icross^ections of 3-body binary formation 
(see lO'Learv et.all 120061 and references therein) . Further- 
more, they included both binary-single and binary-binary 
encounters. Considering the sub-set of BH binaries formed 
in the BH clusters that merge within a Hubble time, they 
determined typically a few AdLIGO detections per year for 
old GCs. 

While such results are remarkable and promising, a 
more detailed study of the dynamics of stellar-mass BHs 
in star clusters with a realistic number of stars is essen- 
tial. As the cross-sections of different processes governing 
the dynamics of BH binaries, viz., 3-body binary formation, 
binary-single star encounters and ionization have different 
dependencies on the number of stars N of the cluster, an ex- 
trapolation to much different N can be problematic. Hence 
it is important to consider clusters with values of N appro- 
priate for massive clusters or GCs. Also, exact treatments 
of the various dynamical processes are crucial for realistic 
predictions of BH-BH merger rates. Finally, a more care- 
ful study of the different dynamical processes leading to the 
formation and evolution of BH binaries in star clusters is 
needed to understand better under which conditions tight 
inspiralling BH binaries can be formed dynamically from a 
star cluster. 

In the present work, we make a detailed study of the 
dynamics of BH-BH binaries formed in a BH sub-cluster, 
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Figure 1. Typical evolution of virial coefficient Q of a star cluster 
with a reflective sphere. Shown is the evolution of the run R3K180 
of Table [T] A rapid initial heating is observed, caused by super- 
clastic encounters between BHs, followed by a saturation of the 
heating curve, caused by the enhanced escape rate of stars due to 
the heating (see text). 



as introduced above. In particular, we investigate whether 
hard enough BH binaries that can merge via gravitational 
radiation in a Hubble time within the cluster or after getting 
ejected from the cluster, can be formed in such a sub-cluster. 
To that end, we perform direct N-body integrations of con- 
centrated star clusters (half-mass radius 1 pc) consist- 
ing of iV sj 10 5 low-mass stars in which a certain number of 
stellar-mass BHs is added, representing a star cluster with 
an evolved stellar population. 

The present paper is organized as follows: In Sec. [2] we 
describe our simulations in detail. We discuss the various ele- 
ments and assumptions of the simulations and summarize all 
the runs that we perform (Sec. 12.2) . We also discuss the use 
of a reflective boundary to simulate only the core of a star 
cluster (Sec. 12. 3) ). In Sec. Owe discuss our results with partic- 
ular emphasis on the dynamical BH binaries formed during 
the simulations. We discuss the BH-BH mergers that occur 
within the clusters and also the merger timescales of the BH 
binaries that escape from the clusters (Sec. I3.2J1 and obtain 
the distributions of merger-times for both cases (Sec. 13. 3[) . 
Finally, in Sec. 2] we interpret our results in the context of 
different types of star clusters that are observed and provide 
estimates of BH-BH merger detection rates. 



2 MODELS 

To study the dynamics of BHs in star clusters, we perform 
direct N-body simulations using NBODY6 with star clusters 
in which a certain number of BHs are added initially. The ini- 
tial density distribution follows a Plummer model with half- 
mass radius ru ^ 1 pc consisting of TV ^ 10 5 low-mass main- 
sequence stars in the mass-range 0.5Mq ^ m ^ l.OM©. Ob- 
served half mass radii for GCs and open clusters are usually 
larger, but as we expect the clusters to expand considerably 
due to the heating caused by the encounters in the BH core 
(see Sec. [TJ , we begin with more concentrated clusters. 

Black holes formed through stellar evolution are typ- 
ically within the mass-range 8Mq < Mbh < 12Mq for 



stellar progenitors with solar-like metallicity (|Casares||2007l ; 
iBelczvnski et all 120091 ). The exact form of the BH mass- 
function is still debated and depends on the metallicity of 
the par ent stellar population a nd the stellar wind mass-loss 
model l|Belczvnski et all l2009h . In a dense stellar environ- 
ment it is further affected by the frequent stel lar merg- 
ers and binary coalescence l|Belczvnski et al1l2008l ). In view 
of the uncertainty of the BH mass-function, we consider 
only equal-mass BHs, with Mbh = WMq in this work. 
Such a value was also assumed by earlier authors, e .g. , 
iPortegies Zwart fc McMillan! (|200fj| ). lBenacciuistal |2002ri . 

A specified number of BHs are added to each cluster, 
and we assume that the BHs follow the same spatial distribu- 
tion as the stars initially. The number of BHs added depends 
on the BH retention fraction. We explore both full retention 
and the case where half of the BHs are ejected from the clus- 
ter by natal kicks. With such a cluster, we mimic the epoch 
at which massive stars have already evolved and produce 
BHs. While such a cluster is not completely representative 
of the more realistic case in which the BHs are produced 
from stellar evolution during the very early phases, it still 
serves the purpose of studying the dynamics of the BHs since 
this becomes important only later after segregation of the 
BHs into the cluster core (see Sec. [3| . 

We perfo rm our simulations using NBODY6 code 
l|AarsethlB003l ) enabled for use with Graphical Processing 
Units (GPU). We have made arrangements in NBODY6 to 
put a specified number Nbh of BHs of a given mass Mbh by 
picking stars randomly from the cluster and replacing them 
by the BHs. Necessary rescalings have been done to compen- 
sate for the excess mass gained by the cluster. These runs are 
performed on NVIDIA 9800 GX2/GTX 280/GTX 285 GPU 
platforms, located at the Argelander-Institut fur Astronomic 
(AlfA), University of Bonn, Germany, with the GPU en- 
abled version of NBODY6 (seel verre Aarseth's Homepage : 
Isverre.comh . 



2.1 BH-BH mergers 

To evolve the BH-BH bin aries due to GW emission, Pe- 
ters' formula (|Peterslll964r ) is utilized in NBODY6 (also see 
iBaumgardt et all 20061 ) .which provides approximate orbit- 
averaged rates of evolution of binary semi-major axis a and 
eccentricity e due to GW emission. According to this for- 
mula, the merger time T mrg of an equal-mass BH-BH binary 
due to GW emission is given by 



T m ^I50Myr(^_) £ 



(1-e 



2n7/2 



(1) 



Peters' formula is limited up to the mass quadrapole 
terms of the radiating system. Close to the merger, higher 
order PN terms become im portant, which m odifies the above 
value of the merger time ([Blanchct 2006). In the present 
work however, we are primarily interested in an overall 
statistics of the merger events rather than the detailed or- 
bit of a single BH-BH merger, so that the higher order PN 
corrections are not crucial. Hence we restrict ourselves to Pe- 
ters' formula. Moreover, the orbit of a tight BH-BH binary 
gets modified by dynamical encounters, which will anyway 
modify its merger time to a much larger extent compared to 
that for the unperturbed binary with the higher PN terms. 
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In NB0DY6, the orbital evolution of compact binaries 
is also considered when the binary is inside a hierarchy. Thus 
a tight BH binary will continue to shrink even if it acquires 
an outer member forming a hierarchical triple, which can 
often happen due to the strong focussing effect of BH bi- 
naries. Also, energy removal due to GW (bursts) during a 
close hyperbolic passage between two BHs is considered in 

the code. 

Numerical simulations of BH-BH mergers (see iHughej 

2009 for an excellent review) indicate that for unequal-mass 
BHs or even for equal-mass BHs with unequal spins, the 
merged BH product acquires a velocity-kick of typically 100 
km s _1 or more due to an asymmetry in momentum out- 
flow from the system, associated with the GW emission. 
Although we consider equal-mass BHs, the merger kicks as- 
sociated with the inequality of the spins of the merging BHs 
would be generally sufficient to eject merged BHs from the 
cluster. Therefore, in our simulations we provide an arbi- 
trarily large kick of 150 km s _1 immediately after a BH-BH 
merger, to make sure that the merged BH escapes. 

2.2 Computations 

To study the rate of BH-BH mergers coming from a star clus- 
ter, we perform simulations of isolated star clusters with sin- 
gle low-mass stars and BHs as mentioned above. The forma- 
tion of the BH-core through mass segregation and its dynam- 
ics remain largely unaffected by the presence of a tidal field, 
which mainly affects stars near the tidal boundary. While 
the enhanced remov al of stars acce lerates the core-collapse of 
the cluster (see e. a. . ISpitzerlll987l . Ch. 3), the latter is more 
strongly enhanced by the collapse of the BHs themselves 
(~ 100 Myr timescale, see below), so that the effect of any 
tidal field is only second-order. Hence, isolated clusters are 
good enough for our purposes. Further, for simplicity, we do 
not take into account primordial binaries in this initial study. 
The primordial binary fraction in GCs and their period 
distribution is still widely debated llAshman fc ZepJ fl998l ; 
Bcllazzini et al ] |2002l :ls ommariva et aDT 2009). The presence 
of a primordial binary population can however significantly 
influence the dynamics of stellar-mass BHs which we defer 
for a future more detailed study. 

For sol ar-like metallicity, Egglet o n's s tellar evolu- 
tion model (|Eggleton. Fitchett fc Toutl (|l989h . adapted in 
NBODY6) gives about N BH ~ 200 BH s for a cluster with 
TV = 10 s stars following a Kroupa IMF (lKroupall200lT l. The 
above Nbh (or its proportion with TV) is thus an upper limit 
to the number of BHs in a GC that corresponds to a full re- 
tention (i.e., no or low natal kicks for all BHs). 

We perform 2 simulations with TV = 4.5 x 10 4 and 
Nbh = 80, 2 runs with TV = 6.5 x 10 4 , N BH = 110, i.e., 
about full BH-retention. Two of the above runs are repeated 
with half the above Nbhs. Also, one run with TV = 5 x 10 4 
and excess Nbh = 200, appropriate for a top-heavy IMF has 
been performed. Finally, we do 2 runs with N = 10 5 with 
Nbh = 80 (about 50% retention fraction) and 200 (full re- 
tention). All the clusters consist of low- mass stars between 
O.5M < m < l.OM with a Kroupa IMF and the BHs have 
Mbh = 10M© , as discussed in detail in the beginning of the 
section. 

In addition to these systems, we also study stellar-mass 
BHs in clusters with smaller TV representing open clusters, 
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Figure 2. Typical example of the mass segregation of BHs. 
Shown is the radial position R vs. time t (top panel) of all BHs for 
model C50K80 of Table [T] Each of the points represents a BH. 
The BHs segregate within 50 Myr during which Nbh remains 
unchanged (bottom panel). As the BH sub-cluster becomes dense 
enough that BH-BH binary formation takes place through 3-body 
encounters, BHs and BH binaries are ejected from the BH-core 
and Nbh starts to decrease. 



in order to estimate the lower-limit in cluster mass for the 
occurrence of BH-BH mergers. We perform 10 runs with 
N = 5 x 10 3 , N = 10 4 and N = 2.5 x 10 4 each with full BH- 
retention (i.e., Nbh = 12, 20 and 50 respectively). Results 
of all our runs are summarized in Table [T] 



2.3 Simulation of a GC core: reflective boundary 

We also perform simulations with a smaller number of stars 
and BHs that are confined within a reflecting spherical 
boundary. With such a dynamical system one can mimic 
the core of a massive cluster, where the BHs are concen- 
trated after mass segregation. The advantage of this ap- 
proach is that one can simulate the evolution of a massive 
cluster with much fewer stars. We simulate TV = 3000 — 4000 
stars packed within 0.4 pc. This provides a stellar density 
of ~ 10 4 Mq pc -3 , appropriate for the core-density of a 
massive cluster. The initial BH population is taken to be 
Nbh ~ 100 or 200, representing half or full BH-retention 
respectively, of a TV = 10 5 star cluster. In this way, the sim- 
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ulation of the BH-core of a massive GC can be performed 
much faster than the simulation of the whole cluster. 

To implement the reflective boundary, we simply con- 
sider all outgoing stars beyond the reflective sphere R a 
within each block and reflect them elastically. This is done 
in the main integratio n loop. For th ese stars, we compute 
the force polynomials ljAarsethll2003l ) separately (using the 
CPU) which makes the code slower by factors of 2-3 depend- 
ing on the stellar density. While the behavior of a N-body 
code in presence of a reflective boundary needs to be stud- 
ied in more detail, the present implementation seems to be 
fairly stable with single-stars and BHs, with energy errors 
marginally larger than that for a free cluster. (In fact, in the 
NBODY6 package, some routines are already available for 
optionally implementing the reflective boundary). For stars 
and BHs beyond a pre-assigned speed v esc , representing the 
escape speed of the host cluster, we do allow them to escape 
through the reflective boundary. 

It is of course important to note that a cluster sim- 
ulation within a reflective boundary cannot represent the 
evolution of a free cluster, in particular, since the cluster ex- 
pands due to the BH-heating which continually reduces the 
stellar density in the core. With the reflective boundary, the 
stellar density still decreases because of the escape of stars 
with speed v > v esc , but it generally does not mimic the 
cluster expansion. The cooling effect of these escaping stars 
also helps to inhibit the runaway heating of the cluster due 
to the binding energy released by binary-star/binary-binary 
encounters. The escaping single BHs and BH binaries are 
particularly efficient in this mass-loss cooling due to their 
larger masses. Fig. [T] shows a typical example from one of 
our runs (R3KI80 of Table [1] see below) in which the be- 
havior of the virial coefficient Q = T/V (T = total kinetic 
energy of the system, V = total potential energy) is shown. 
While Q initially increases quickly indicating rapid heating 
of the cluster, the escape rate of stars and BHs also increases 
due to the increased velocity dispersion. The latter effect in 
turn increases the cooling rate so that Q finally tends to 
flatten to a value for which the dynamical heating is bal- 
anced by the escape-cooling of the cluster. For the reflective 
boundary clusters used in our simulations (see above), the 
initial velocity-dispersion of the stars within the reflective 
boundary is chosen to be vo fa 7 km s" 1 and escape speed 
v esc ~ 24 km s _1 . Typically, during each run, the virial co- 
efficient Q increases by a factor of about 6.5 and most of the 
growth takes place during the first few hundred Myr, as in 
Fig. [T] As v ~ Q 1//2 , it increases to w 17 km s _1 , which is 
typical for the central velocity dispersion of a massive clus- 
ter. 

In a reflective boundary cluster, the BH binaries start 
forming immediately after the start of the integration since 
the system already initiates with a BH-core. We perform a 
set of 6 runs with reflective boundary clusters (see Table [TJ. 
We shall discuss results of our reflective boundary simula- 
tions in Sec. [3] 

3 RESULTS: BLACK HOLE BINARY 
MERGERS AND ESCAPERS 

In this section we discuss the results of the simulations in- 
troduced in Sec. 12.21 As already mentioned, we focus on 



tight BH-BH binaries both within the cluster and the es- 
caped ones that can merge within a few Gyr. Table [T] pro- 
vides an overview of the results of our simulations for both 
isolated and reflective clusters. It shows the number of BH- 
BH mergers within the cluster for each computation and 
the corresponding merger times. Also, for each computation, 
the numbers of escaped BH-BH binaries that merge within 
3 Gyr are shown. For convenience of discussion, we utilize 
one of these models, viz., the one with identity C50K80 (see 
Table [TJ — all others generally posses similar properties. 

Fig. [2] (top panel) demonstrates the mass segregation 
of the BHs in the cluster which takes about 50 Myr. Before 
the BHs segregate to within about 0.3 pc of the cluster cen- 
ter, the BH-density is small so dynamical encounters among 
the BHs are not significant and Nbh remains constant. As 
the BHs segregate within about 0.3 pc of the cluster core, 
the BH-density of this sub-cluster becomes high enough to 
form BH-BH binaries through 3-body encounters. Once BH- 
BH binaries start forming, single and BH binaries begin to 
escape from the BH-core by the recoils received due to re- 
peated super-elastic encounters between the single and bi- 
nary BHs (see Sec.[T|. In Fig. [2] (top panel) one can clearly 
distinguish the two phases of the BH subsystem — the ini- 
tial segregation phase and the BH-core formation, the radial 
positions of the BHs in the latter phase being scattered out- 
wards due to the recoils. The decrease of Nbh during this 
phase is also shown in Fig. [5] (bottom panel) . 

The time i SC g taken by the BHs to segregate to the 
cluster core to form the BH sub-cluster is essentially the 
core-coll apse tim e for t he initial BH cluster. t seg is given by 
(see e.q.. lSpitzed (|l987f ), Ch. 3): 

where, t cc is the core-collapse time of the host star cluster it- 
self (i.e., without the BHs) and (m) is its mean stellar mass. 
According to numerical experimen ts by several authors {e.g., 
iBaumgardt, Hut fe Heggie] l|2002h 'l. a Plummer cluster takes 
about t cc « 15i r h (0) to reach the core-collapse stage, where 
trh(0) is the initial half-mass relaxation time of the Plum- 
mer cluster. F or the a b ove e xample, t r h{0) ~ 73 Myr (see 
Eqn. (2-63) of ISpitzerl ljl987t i) and (ra) « 0.6 which gives 
tseg ~ 65 Myr from Eqn. @. This roughly agrees with the 
formation time of the central dense BH-cluster as observed 
in the N-body integration (see Fig. 0. 

3.1 Dynamical black hole-black hole binaries 

A particular dynamically formed BH-BH binary can initially 
be very wide, sometimes with semi-major-axis a more than 
5OOO7?0. However, it then shrinks very rapidly due to re- 
peated super-elastic encounters with the surrounding BHs. 
The collisional hardenin g rate is given by he oc a 2 (see 
iBaneriee fc Gh osh 2006 and references therein). The vari- 
ation in a and e is random due to the stochastic nature 
of t he encounters, but t he binary hardens on average (see 
e.o. jHeggie fc Hutll2003l . Ch. 19 & Ch. 21) and the overall 
hardening rate decreases with a. These encounters include 
both flybys and exchanges with single BHs. When the binary 
is wide, and the recoils it receives are low, so that it usu- 
ally remains in the cluster. However, as it becomes harder, 
the recoils become stronger, and the binary is ejected from 
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Figure 3. Semi-major-axis a vs. evolution time t for all BH-BH pairs formed for the model C50K80, where a different symbol is used 
for each newly formed BH-BH binary. The sequences of the points indicate that each BH binary initially hardens quickly to smaller 
values of a where the hardening rate decreases. The estimated upper-limit a eBC for which a BH-BH pair can escape due to recoil from 
an encounter with a single BH is indicated by the horizontal line. Hardening rates of all BH pairs terminate below this line. Note that 
different traces of the same symbol represent different hardening sequences following encounters. 



the BH-core, but returns to the BH-core due to dynami- 
cal friction (see Sec. [TJ. Finally, the recoil becomes strong 
enough that the B H binary escapes from the cluster (see also 
lBenacauistall2002h . 

For equal mass binary-single star encounters, 40% 
of the binding energy of t he binary is released per 
close encounter on average (ISpitzerl Il987l . Eqn. (6-25); 
iHut. McMillan fc Roman! Il992h . This leads to the follow- 
ing expression for the average recoil velocity of a BH binary 
due to encounters with single BHs: 

(Vrec) = 6.75 x 1(T 2 ° MbH . (3) 
a 

Hence, the BH binary can be ejected due to recoil from a 
cluster with escape velocity v esc for a < a esc , where a esc is 
obtained by setting the left hand side of Eqn. (|3]) to v\ ac . 
For the above Plummer clu ster, we then get from the value 
of its central potential (see iHeggie fc Hutl [20031 . Table 8.1) 
a eac m 500Rq , which is shown as solid line in Fig. [3] The 
escaped binaries in our calculations are generally found to 
be harder than a esc except for a few systems, which escape 
within triple-BHs. 

Fig. [3] provides an overall impression of the hardening 
of the BH binaries where the values of a for all BH-BH 
binaries formed during the computation are plotted. Each 
newly formed BH pair is represented by a different symbol. 
From the sequences of points, it can be seen that each of the 
BH-BH pairs exhibits an overall tendency of rapid hardening 
and the steep a(t) dependence tends to flatten as a decreases. 
The termination of a particular curve generally indicates 
escape of the corresponding BH binary from the cluster (or 
more rarely a merger, see Sec. 13 . 2 p . Typically, a BH binary 
leaves the cluster with a < 2007?© (see Fig. [5}. 



3.2 Mergers and escapers 

Can a BH-pair be hardened to small enough a (and eccentric 
enough at the same time) for it to merge via GW radiation? 
To investigate this question, we consider the positions of the 
BH-BH binaries within the cluster in a a vs. (1 — e 2 ) plane 
(Fig. [3), where each newly formed BH pair is represented 
by a different symbol as in Fig. [3] The points are color- 
coded with the evolution time in Myr, the color-scale being 
displayed on the right. Although for each BH-pair a and e 
fluctuate over the plane, the changes occur on a collision 
time-scale which is ~ Myr. Since the orbital period of the 
binaries corresponding to these points is much shorter (from 
~ days to years), these points generally represent binaries 
which are stable over many orbits. Overplotted on Fig. [3J 
are lines of constant GW merger time, T mrg , as given by 
Eqn. {TJ. While most of the points have very large merger 
time, a few of them do lie around the T mrg — 10 Myr line. 
This indicates that these binaries are indeed hardened up to 
small enough a and/or acquire sufficient eccentricity that if 
they are left unperturbed, they can merge via GW emission 
within several Myr. This feature is found to be generally true 
for all the computations reported here (see Table [TJ, ie, each 
of them does produce BH-pairs that are capable of merging 
within several Myr. It is important to note that the short 
merger timescale of these binaries is due to their high eccen- 
tricities since they are generally far too wide (~ WORq) to 
merge unless they are highly eccentric (see Eqn. Q). This 
is also true for the escaping BH binaries (see below). 

However, these merging BH-pairs, particularly due to 
their large eccentricity and focussing effect, can still be per- 
turbed by further encounters on timescales ~ Myr which 
can often prevent them from merging. In the particular ex- 
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Figure 4. Positions of all the BH binaries within the cluster in a 1 — e 2 vs. a plane where different symbols are used to distinguish 
between different BH pairs (for C50K80). The color-coding of the points, as indicated by the color-scale, represents the evolution time 
(in Myr) of the cluster at which they appear at a particular location in the above plane. See text for details. 



ample shown in Fig. [4] only one of the "merging candi- 
dates" (Tmrg 10 Myr) could actually merge within the 
cluster (inverted triangles) . For all the other candidates, in- 
cluding the one which escaped (filled dots), the eccentricity 
became much smaller due to encounters so that the merger 
timescales increased considerably. 

As for the escapers, since they remain unperturbed af- 
terwards, all with GW merger times smaller than a Hub- 
ble time are of interest. However, we find that typically for 
each cluster there are relatively few BH binaries with GW 
merger times T mrg > 3 Gyr among those that merge within 
a Hubble time (see e.g., Fig. [5] and Fig. [6}. In the above 
example, i.e., model C50K80, there is one BH-BH escaper 
with T m r a ~ 1 Gyr. A more interesting example is shown 
in Fig. [5] corresponding to the model C65K110 of Table [1] 
where 2 binaries merge in m 3 Gyr, one in ~ 1 Gyr and one 
in < 10 Myr. 

The results of the present computations, as summarized 
in Table [JJ indicate that a medium-mass to massive cluster 
with sufficient BH-retention is likely to have at least one BH- 
pair that can merge within the cluster within a time-span of 
a few Gyr. Each of these clusters also typically eject a few 
BH binaries that can merge within a Hubble time. In one of 
the models (C50K200), one BH-BH hyperbolic collision also 
occurred at t mrg ~ 820 Myr. As such hyperbolic mergers 
did not show up in any of the other models, they must be 
much rarer than GW driven mergers. 



3.3 Merger-time distribution 

As a representation of the merger time distribution of merg- 
ers within the cluster, a combined distribution of t mrg of all 
the BH binaries that merged within the clusters is shown in 
Fig. [6] (top panel). This is justified, as the values of BH-BH 



binary merger times t mrg within the cluster, shown in Ta- 
ble!]] apparently do not indicate any correlation of t mrg with 
the cluster parameters, in particular the number of stars N . 
Both the isolated clusters and the reflective boundary clus- 
ters are included. The t mrg distribution in Fig. [6] indicates 
that Nmr g decreases with t mrg . This can be expected, as at 
later times, Nbh in a cluster decreases, so the hardening 
rate of BH binaries also decreases, making their merger less 
probable. Comparing our results for clusters with different 
values of N , we find that the rate of depletion of BHs from 
the BH core during the first few Gyr does not depend appre- 
ciably on N. In our computations, BHs are depleted within 
a few Gyr for most of the models. Clusters with large Nbhs, 
e.g., C100K200, C50K200 do indicate longer BH-retention, 
but this is due to the expansion of the cluster (and therefore 
the BH sub-cluster itself) because of the heating by the BH- 
core (see Sec. [lj and hence it is unlikely that tight BH-BH 
binaries can be formed during the later phases. 

Fig. [5] (bottom panel) shows the distribution of the 
merger times for the escaping BH binaries from all the com- 
putations of Table [JJ Note that the number of tight esca- 
pers in Table [JJ also does not indicate any dependence on 
N. For each escaper f mr9 includes its time of escape t esc , 
i.e., t mrg = t eac + T mrg , where T mrg is calculated from 
Eqn. ([JJ) using the values of a and e with which it escaped. 
The merging rate among escaped binaries also shows a de- 
crease with time. Such a decrease is expected from Eqn. {JJ. 
For example, if we consider that the escapers that merge 
within a Hubble time have a mean radius a ~ 50Rq and a 
thermal eccentricity distribution dN mrg /de oc 2e, it is easy 
to show from Eqn. (UJ that dN mT g/dT mrg oc Tmr g 7 . The 
above merger-time distribution indicate that significantly 
more BH-BH mergers occur outside the clusters, i.e., among 
the escapers, than within the clusters. We shall return to this 
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Table 1. Summary of the computations performed for isolated clusters and those with reflective boundary, see Sec. I2.2I &; Sec. 12.31 The 
meaning of different columns is as follows: Col. (1): Identity of the particular model — similar values with different names (ending with 
A,B etc) imply computations repeated with different random seeds, Col. (2): Total number of stars N, Col. (3): Number of simulations 
Nsim with the particular cluster, Col. (4): Initial half- mass radius of the cluster r h (0) (isolated cluster) or radius of reflective sphere 
R s , Col. (5): Initial number of BHs Ngg(0), Col. (6): Total number of BH-BH binary mergers within the cluster N mrg , Col. (7): The 
times t mrg corresponding to the mergers, Col. (8): Number of escaped BH-pairs iV esc — the three values of N eac are those with T mrg < 
3 Gyr, 1 Gyr and 100 Myr respectively, Col. (9): BH-BH Merger rate 7?-AdLlGO detected by AdLIGO assuming that the corresponding 
model cluster has a space density of p c ; = 3.5 h 3 Mpc -3 (see Sec. 14. ill . 



Model name 


N 




r h (0) or R s (pc) 


N BH (0) 


Nmrg 


tmrg (Myr) 


Nesc 




Isolated clusters 


C5K12 


5000 


10 


1.0 


12 











C10K20 


10000 


10 


1.0 


20 











C25K50 


25000 


10 


1.0 


50 







3 11 




C50K80 


45000 


1 


1.0 


80 


1 


698.3 


3 10 


28(±14) 


C50K80.1 


45000 


1 


0.5 


SO 


2 


217.1, 236.6 


3 2 1 


35(±15) 


C50K40.1 


45000 


1 


0.5 


40 







111 


7(±7) 


C50K200 


50000 


1 


1.0 


200 


2 


100.8, 467.8 





14(±10) 


C65K110 


65000 


1 


1.0 


110 


1 


314.6 


4 2 1 


35(±15) 


C65K110.1 


65000 


1 


0.5 


110 







4 3 1 


28(±14) 


C65K55.1 


65000 


1 


0.5 


55 


1 


160.5 


10 


14(±10) 


C100K80 


100000 


1 


1.0 


80 


2 


219.4, 603.2 


5 2 1 


42(±15) 


C100K200 


100000 


1 


1.0 


200 







5 4 4 


28(±14) 


Reflective boundary 


R3K180 


3000 


1 


0.4 


180 


1 


1723.9 


5 3 1 


35(±15) 


R4K180A 


4000 


1 


0.4 


180 


1 


3008.8 


2 2 1 


21(±12) 


R4K180B 


4000 


1 


0.4 


180 


2 


100.2, 1966.5 


2 10 


28(±14) 


R3K100 


3000 


1 


0.4 


100 


2 


3052.8, 3645.9 


110 


18(±10) 


R4K100A 


4000 


1 


0.4 


100 


2 


104.4, 814.2 


3 3 1 


28(±14) 


R4K100B 


4000 


1 


0.4 


100 


1 


1135.3 


3 3 3 


28(±14) 



CD 



0.01 



0.001 




300 



a(Fy 



Figure 5. Positions of all escaped BH binaries in a 1 — e 2 vs. a 
plane for the model C65K110. Lines of constant merger times are 
plotted as in Fig. [4] 



4 DISCUSSION 

Our present study indicates that centrally concentrated star 
clusters, with > 4.5 x 10 4 are capable of dynamically 
producing BH binaries that can merge within a few Gyr, 
provided a significant number of BHs are retained in the 
clusters after their birth. The results of our simulations (see 
Table[TJ imply that most of the BH-BH mergers occur within 
the first few Gyr of cluster evolution for both mergers within 
the cluster and mergers of escaped BH binaries. 

The above results imply that an important class of can- 
didates for dynamically forming BH binaries that merge 
at the present epoch are star clusters with initial mass 
M ct > 3 x 10 4 A/(rQ, which are less than few Gyr old. 
Such clusters represent intermediate-age massive clusters 
(hereafter IMC) with initial masses close to the upper- 
limit of the initial cluster mass function (ICMF ) in spiral 
JWiedner. Kroupa fe Larsenll2004l ; lLarsen|[2009a| ) and star- 
burst ( Gieles et al.ll2006l ) galaxies. While it is not impossi- 
ble to obtain BH-BH mergers within a Hubble time from 
lower-mass clusters, the overall BH-BH merger and escape 



point in Sec. [4] The median merger time for BH-mergers 
within the clusters is t mrg » 1 Gyr and for the escaped BH 
binaries it is ~ 3 Gyr. 



1 To correlate with the observed clusters, we consider the mass 
M c i of the parent cluster, i.e., the cluster mass before the BHs 
are formed through stellar evolution, which are heavier than the 
clusters of low mass stars that we model, for the same value of 
N . For a given simulated cluster, we estimate the corresponding 
parent mass M c i by weighting the mean stellar mass of (m) c i rj 
O.6M of a Kroupa IMF with star from O.IMq upto 100Af Q with 
the total number of stars N for that cluster. 



Stellar-mass black holes in star clusters 9 




2000 4000 6000 8000 10000 12000 14000 

Wg ( M y) 



14 
12 
10 
8 
6 
4 
2 




- 1 T— 



_l L 



m rm 



2000 4000 6000 8000 10000 12000 14000 
Wg (Myr) 

Figure 6. Top: Distribution of the merger times t mrg for BH 
binary mergers within the cluster for the models of Table HI Bot- 
tom: Distribution of the merger times tmrg for escaped BH bi- 
naries for the models of Table [T] (see text). 



rates strongly decrease with cluster-mass as Table \T\ indi- 
cates. For M ci < 1.5 x 1O 4 M , mergers already become 
much rarer (see Table [TJ. Due to the statistical nature of 
merger or ejection events it is ambiguous to set any well- 
defined limit on the cluster-mass beyond which these events 
become appreciable (such an estimate would also require a 
much larger number of N-body integrations). In view of our 
results, Md ~ 3 x 1O 4 M0 is a representative lower limit be- 
yond which an appreciable number of mergers and escapers 
merging within a Hubble time can be obtained. 

Old globular clusters, which can be about 10 times or 
more massive, are expected to produce mergers or escapers 
more efficiently. As the timescale of depletion of BHs from 
the BH cluster is nearly independent of the parent cluster 
mass (see Sec. I3.3|l . GCs can also be expected to produce 
BH-BH mergers over similar time-span as the IMCs, i.e., 
within the first few Gyr of evolution. Since GCs are typ- 
ically much older (~ 10 Gyr), they do not contribute sig- 
nificantly to the present-day merger rate, since most of the 
mergers from them would have occurred earlier. Consider- 
ing the light-travel time of « 4.5 Gyr from the maximum 
distance D w 1500 Mpc form which these BH-BH binaries 
can be detected by "AdLIGO"(see below), only GCs close 
to the above distance could contribute detectable events, 
mostly from escaped BH-BH binaries. 



On the other hand, young massive clusters with ages less 
than 50 Myr, repr esenting star c lusters near the high-mass 
end of the ICMF |l arsenl [2009b! ') . are generally too young 
to produce BH-BH mergers. All models in Table [1] produce 
mergers significantly later than this age (except one escaped 
BH binary in each of the models C65K110 and C100K200). 
Hence, IMCs seem to be most likely candidates for producing 
observable BH-BH mergers dynamically. 



4.1 Detection rate 

We now make an estimate of the BH-BH merger detec- 
tion rate from IMCs by ground-based GW observatories 
like LIGO and AdLIGO. In estimating the overall BH-BH 
merger rate using the results of our model clusters, one needs 
to consider the distribution of the cluster parameters that 
are varied over the models, viz., cluster mass, half-mass ra- 
dius, and BH retention fraction. Such distributions are far 
from being well determined, except for the mass di stribution 
for young clusters in spiral and starburst ga laxies l|Bik et al.l 
2003; Bas tian fc Lam crs 20031; lGielesi r2004). Therefore, de- 
termination of an overall merger rate considering the distri- 
bution of our computed clusters can be ambiguous. Hence, 
as a useful alternative, we determine the BH binary merger 
detection rates for each of the cluster models in Table[T]that 
gives an appreciable number of mergers, for a representative 
density of IMCs. Su ch an approach has been considered by 
earlier authors, e.q.. lO'Learv et.all (|2006l ) and can provide a 
reasonable idea of the rate of detection of BH-BH mergers 
from IMCs. 

As an estimate of the space density of IMCs, 
we adopt that for yo u ng p opulous clusters in 
IPortegies Zwart fc McMillan! (120001 ) ■ which has a simi- 



lar mass-range as the IMCs: 



Pel 



3.5 hr Mpc" 



(4) 



where h is the Hubble parameter, defined as Hp / 100 km s 1 , 
Ho being the Hubble constant <|PeebIedll993h . The above 
space density has been derived from the spac e densities 
of sp iral, blue elliptical and starburst galaxies l|Hevl et al.l 
Il997t) assuming that young populous c lusters have the same 
specific frequencie s (Sn) as old GCs (|van den Berghlll995l ; 
iMcLaughlinll 19991 ). but in absence of any firm determination 
of the Sns of the former. We compute the detection rate for 
each model cluster assuming that it has a space density of 
the above value. 

The LIGO/ AdLIGO detection rate of BH-BH merg- 
ers from a particular m odel cluster can be estimated from 
(jBelczynski et al.ll2007l and references therein) 

4 3 

TZlIGO = g 71 "^ PclUmrg, (5) 

where lZ mrg is the compact binary merger rate from a cluster 
and D is the maximum distance from which the emitted GW 
from a compact-binary inspiral can be detected. D is given 
by 



D = D 



5/6 



(6) 



where D = 18.4 and 300 Mpc for LIGO and AdLIGO re- 
spectively. The quantity M c h is the "chirp mass" of the com- 
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pact binary with component masses mi and mz, which is 
given by 



M rh , = 



(m 1 m 2 ) 3/5 
'mi + m 2 ) 1 / 5 ' 



(7) 



and M c h,nsns = 1.2A4q is that for a binary with two 1.4Mq 
neutron stars. 

For a BH binary with mi = m 2 = IOMq, M c h = 
8.71M Q which gives D « 1500 Mpc for AdLIGO. The 
AdLIGO detection rates 7^AdLiGO (mean over 3 Gyr taking 
into account the time of escape t esc of the escaped binaries) 
for the model clusters are shown in Table [T] where the cur- 
rently accepted value of the Hubble parameter h = 0.73 is 
assumed. The error in each detection rate is simply obtained 
from the Poisson dispersion of the total number of mergers 
for the corresponding cluster. Note that these detection rates 
are for clusters with solar-like metallicity which is implied 
by our assumption of WMq BHs for all the clusters. 

To obtain a basic estimate of the overall detection rate 
of BH-BH mergers from IMCs, we consider the subset of 
our computed models that are isolated clusters with full 
BH retention (see Table [1} . We take the mass function of 
the IMCs to be a power law with index a = —2 which 
is the approximate in dex of the ICMF in spira l and star- 
burst galaxies (see e. a. . lGielesll2004l ; lLarsenll2009al ') . Then the 
weighted average of the corresponding AdLIGO detection 
rates is 7?.AdLiGO ~ 31(±7) yr _1 , which estimates the total 
present-day detection rate of BH-BH mergers from IMCs ex- 
pected for AdLIGO. The corresponding LIGO detection rate 
is negligible, 7^ligo w 7.4 x 10" 3 yr" 1 . Note that these BH- 
BH detection rates are only lower limits. First, the observed 
population of star clusters can be an underestimation by a 
factor of 2 owing to their dissolution in the tidal field of their 
host galaxies fsee lPortegies Zwart fc McMillanll2000] and ref- 
erences therein). Second, the above detection rates are only 
from IMCs and there can be additional contributions from 
GCs (see above). 

In comparing the AdLIGO detection rates from our 
computations with those from earlier works, we note that 
our rates are typically an order of mag n itude smaller than 
those of iPortegies Zwart fc McMillan! (|200Ch . but about 
one order of magn itude larger than those obtained by 
lO'Learv et.a 3(120061). The principal origin of the former dif- 
ference is due to the fact that while we considered only 
IMCs distinguishing them as the most appropriate can- 
didates for producing present-day B H-BH mergers (see 
above), IPortegies Zwart fc McMillan! (|2000h also included 
GCs, which have considerably larger spatial density and also 
larger fraction of BH-BH binaries merging within the Hubble 
time, as obtained from their analytic extrapolations. Note 
that the number of escapers as obtained by them for young 
populous clusters and the fra ction of them merging within 
a Hu bble time (see Table 1 of IPortegies Zwart fc McMillan! 
2000) is similar to that obtained from the present compu- 
tations, implying qualitative agreement. The above authors 
apparently did not consider the time scales of formation and 
depletion of the BH sub-syste ms in their pre l imina ry study. 
On the other hand, although lO'Learv et.all (|2006l ) consid- 
ered clusters significantly more massive than our's in their 
Monte-Carlo approach, so that larger merger detection rates 
can be expected, their clusters were much older (8 Gyr and 



13 Gyr) than IMCs which, in accordance with our results, 
accounts for their much lower detection rate. 

It is interesting to note that the dynamical BH-BH 
merger detection rates obtained by us are typically an order 
of magnitude higher than that from primord ial stellar bina- 
ries as predicted by iBelczvnski et~a l. (2007) based on their 
revised binary evolution model, and is similar to that for the 
isolated NS-NS binaries derived by them. Hence our results 
imply that dynamical BH-BH binaries constitute the domi- 
nant contribution to the BH-BH merger detection. Thus, the 
dynamical BH-BH inspirals from star clusters seem to be a 
promising channel for GW detection by the future AdLIGO, 
although their estimated detection rate with the present 
LIGO detector is negligible, in agreement with the hitherto 
non-detection of GW. 



4.2 Limitations and outlook 

The work presented here is a first step towards a detailed 
study of the dynamics of stellar-mass BHs in star clusters 
and the consequences for GW-driven BH mergers, and im- 
provements in several directions are possible. First, we do 
not consider the initial phase of the cluster here when BHs 
form through stellar evolution and a more consistent ap- 
proach would be to begin with a star cluster with a full stel- 
lar spectrum and produce the BHs from stellar evolution. 
Note however that in the present study we have inserted 
numbers of BHs in our clusters of low-mass stars similar to 
what would have been formed from the stellar evolution of a 
cluster following a Kroupa IMF (see Sec. 12. 2p . Stellar evolu- 
tion also produces NSs (about twice in number as the BHs) 
which also segregate to the central region of the cluster. It 
is interesting to study the dynamics of the NS cluster and 
how it is affected by the (more concentrated) BH cluster — 
in particular the formation of tight inspiralling NS-NS and 
NS-BH binaries, which are important for both GW detection 
and GRBs. 

Another aspect that we do not consider in our present 
study is the effect of primordial stellar binaries. Tight stel- 
lar binaries aid th e formation of compa ct binaries through 
double-exchanges (|Grindlav et.alj 2006). in addition to the 
3-body mechanism (see Sec.[TJ which can increase the num- 
ber BH-BH (also BH-NS and NS-NS) binaries formed and 
hence the merger rate. Thus the study of BH-BH binaries in 
star clusters with primordial binaries is an important next 
step. Such studies are in progress and will be presented in 
future papers. 

Finally, the number of N-body computations in this ini- 
tial study is not enough to obtain the BH-BH merger rate as 
a function of cluster mass and BH retention fraction with a 
reasonable accuracy. There are typically one integration per 
cluster model. To obtain merger rate dependencies with clus- 
ter parameters, e.g., mass, half mass radius, binary fraction 
and BH retention, many computations are needed within 
smaller intervals, which involves much larger time and com- 
puting capacity than that utilized for the present project. 
Such results, combined with improved knowledge of cluster 
parameter distributions and BH formation through super- 
nova explosions of massive stars, would provide a robust 
estimate of the BH-BH merger detection rate. Conversely, 
with such a merger rate function, the detection of BH-BH 
mergers by GW detectors like AdLIGO in the near future 
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will shed light on the above mentioned long-standing ques- 
tions. 



ACKNOWLEDGMENTS 

We are thankful to Sverre Aarseth of the Institute of As- 
tronomy, Cambridge, UK, for many discussions and sugges- 
tions and his timely modifications of the NBODY6 code, 
which have been very helpful for our computations. We 
thank Keigo Nitadori for his development and regular up- 
dates of the GPU libraries utilized in NBODY6. We thank 
the anonymous referee for very helpful comments and sug- 
gestions which have improved several parts of this paper sig- 
nificantly. This work has been supported by the Alexander 
von Humboldt Foundation. 



REFERENCES 

Aarseth, S.J., 2003, "Gravitational N-Body Simulations", 

Cambridge University Press. 
Sverre Aarseth's Homepage: |http://www.sverre.com/| 
Amaro-Seoane, P., Gair, J.R., Freitag, M., Miller, M.C., 

Mandel, L, Cutler, C.J. and Babak, S., 2007, CQG, 24, 

113. 

Ashman, K.M. and Zepf, S.E., 1998, "Globular Cluster Sys- 
tems", Cambridge Univ. Press. 

Banerjee, S. and Ghosh, P., 2006, MNRAS, 373, 1188. 

Bastian, N. and Lamers, H.J.G.L.M, 2003, in "Extragalac- 
tic Globular Cluster Systems", ESO Astrophysics Sym- 
posia, M. Kissler-Patig eds., Springer- Verlag. 

Baumgardt, H., Hut, P. and Heggie, D.C., 2002, MNRAS, 
336, 1069. 

Baumgardt, H., Hopman, C, Potegies Zwart, S. and 

Makino, J., 2006, MNRAS, 372, 467. 
Belczynski, K., Taam, R.E., Kalogera, V., Rasio, F.A., Bu- 

lik, T., ApJ, 2007, 662, 504. 
Belczynski, K., Kalogera, V., Rasio, F., Taam, R., Zezas, 

A., Bulik, T., Maccarone, T., and Ivanova, N. 2008, ApJS, 

174, 223. 

Belczynski, K., Bulik, T., Fryer, C.L., Ruiter, A., Vink, J.S. 

and Hurley, J.R., 2009. larXiv:0904.2784t T. 
Bellazzini, M., Fusi Pecci, F., Messineo, M., Monaco, L. 

and Rood, R.T., 2002, AJ, 123, 1509. 
Benacquista, M.J., 2002, CQG, 19, 1297. 
Bik, A., Lamers, H.J.G.L.M., Bastian, N., Panagia, N. and 

Romaniello, M., 2003, A&A, 397, 473. 
Blanchet, L., 2006, Living Reviews in Relativity, 9, 4. 
Bulik, T., and Belczynski, K., 2003, ApJ, 589, L37. 
Casares, J., 2007, in "Black Holes from Stars to Galaxies 

Across the Range of Masses", V. Karas and G. Matt, eds., 

IAUS, 238, 3. 

Eggleton, P.P., Fitchett, M.J. and Tout, C.A., 1989, ApJ, 
347, 998. 

Gieles, M., 2004, in "The Formation and Evolution of Mas- 
sive Young Star Clusters" , H.J.G.L.M. Lamers, L.J. Smith 
and A. Nota. eds., ASP Conference Series. 

Gieles, M., Larsen, S.S., Scheepmaker, R.A., Bastian, N., 
Haas, M.R. and Lamers, H.J.G.L.M., 2006, A&A, 446, L9. 

Grindlay, J., Portegies Zwart, S. and McMillan, S., 2006, 
Nature Physics, 2, 116. 



Giiltekin, K., Miller, M.C. and Hamilton, D.P., 2004, ApJ, 
616, 221. 

Heggie, D.C., 1975, MNRAS, 173, 729. 

Heggie, D.C. and Hut, P., 2003, "The Gravitational Millon- 
Body Problem: A Multidisciplinary Approach to Star 
Cluster Dynamics", Cambridge University Press, Cam- 
bridge, UK. 

Heyl, J., Colless, M., Ellis, R.S., and Broadhurst, T., 1997, 
MNRAS, 285, 613. 

Hughes, S.A., 2009, to appear in Annual Reviews of As- 
tronomy and Astrophysics, arXiv:0903.4877 

Hut, P. and Verbunt, F., 1983, Nature, 301, 587. 

Hut, P., McMillan, S. and Romani, R.W., 1992, ApJ, 389, 
527. 

Hut, P. et. al, 1992, PASP, 104, 981. 

Kramer, M., 2008, in "Relativistic Astrophysics Legacy and 

Cosmology - Einstein's, ESO Astrophysics Symposia", 

Springer- Verlag Berlin Heidelberg. 
Kroupa, P., 2001, MNRAS, 332, 231. 
Larsen, S0ren. S., 2009a, A&A, 494, 539. 
Larsen, S0ren. S., 2009b, in "Globular Clusters - Guides to 

Galaxies, ESO Astrophysics Symposia", Springer Berlin 

Heidelberg. 

Mackey, A.D., Wilkinson, M.I., Davies, M.B. and Gilmore, 

G.F., 2007, MNRAS, 379, L40. 
McLaughlin, D.E., 1999, AJ, 117, 2398. 
Merritt, D., Piatek, S., Portegies Zwart, S. and Hem- 

sendorf, M., 2004, ApJ, 608, L25. 
OLeary, R. M., Rasio, F. A., Fregeau, J. M., Ivanova, N. 

and OShaughnessy, R., 2006, ApJ, 637, 937. 
O'Shaughnessy, R., Kalogera, V., and Belczynski, K., 2007, 

ApJ, 667, 1048. 
O'Shaughnessy, R., Kim, O, Kalogera, V., and Belczynski, 

K., 2008, ApJ, 672, 479. 
Peebles, P.J.E., 1993, "Principles of Physical Cosmology", 

Princeton University Press, Princeton. 
Peters, P.C., 1964, Phys. Rev. B., 136, 1224. 
Portegies Zwart, S.F. and McMillan, S.L.W., 2000, ApJ, 

528, L17. 

Ransom, S.M., 2008, IAUS, 246, 291. 

Sommariva, V., Piotto, G, Rejkuba, M., Bedin, L.R., Heg- 
gie, D.C, Mathieu, R.D. and Villanova, S., 2009, A&A, 
493, 947. 

Spitzer, L. Jr., 1987, "Dynamical Evolution of Globular 

Clusters", Princeton University Press, 
van den Bergh, S., 1995, AJ, 110, 2700. 
Wiedner, C, Kroupa, P. and Larsen, S.S., 2004, MNRAS, 

350, 1503. 



